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6.7 Thermohaline instability in a porous layer 
-doubly-diffusive instability 
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In Salton Sea geothermal field of Imperial Valley, California, the temperature in the hotwater 
reservoir can be 360°C. Salinity in groundwater is as high as 2.5 x 10 5 ppm. Due to sharp 
contrast of thermal and salt diffusivities (k t ~ 1.5 x 10~ 7 m 2 /s, k$ ~ 1.5 x 10~ 9 m 2 /s), 
convection can result from instability This is called thermohaline convection. 

A very similar mechanism was discovered earlier in oceanography, by Melvin Stern (1960). 
This is another example where essentially the same physics arises in vastly different contexts. 



6.7.1 Heuristic picture 

Monotone instability 
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Consider a fluid lighter at the top and denser at the bottom, with T and C shown in 
Figure 6.7.1. so that dT/dz.dC/dz > 0. 




Figure 6.7.1: Salinity and temperature profiles favoring monotonic instability in a porous 
medium 



Let a fluid parcel rise from z to z + Az. Because the thermal diffusivity K m is much 
greater than the mass difffusivity of salt, k s , the rising water parcel absorbs heat much more 
quickly than it absorbs salt. Hence, it can becomes lighter than the surrounding fluid and 
continue to rise. If a fluid parcel drops by Az, the faster heat loss can make it heavier than 
its sourounding fluid and sinking continues. Thus the fluid system can be unstable. 

Oscillatory instability Consider a cold and fresh water lying over a warm and salty water. 
Assume that the the fluid density decreases with height. The profiles of T and C are shown 
in figure 6.7.2 so that dT/dz, dC/dz < 0. A lighter fluid parcel accidentally moved up by Az 




Figure 6.7.2: Salinity and temperature profiles favoring oscillatory instability in a porous 
medium 



will cool off quickly while losing little salt. Hence, it will come down, resulting in oscillations. 
We shall describe the mathematics of the initial instability theory. 



6.7.2 Governing equations 

Conservation of Fluid Mass: 

V • u = 

Conservation of fiuuid momentum: 

= - Vp - -u + p f g 



Energy: 



Conservation of Salinity : 



Equation of state: 



<9T 

a— + u • VT = K m V 2 T, 
at 

_ (pc) fluid 
(PC) mix 

dC 

n— + u • VC = k s V 2 C 
at 



Pf = p [1 - (3 (T - T ) + (3 S (C - C )} 

Note the sign difference of thermal expansion coefficient and salinity. 
Static equilibrium: 

u = 

T s = r + AT(l-|: 
C s = C + Ac(l-|: 



Ps = Pq~ Po9 



6.7.3 Perturbed flow 

Let u = u' be the velocity disturbance, and let 

t = t s + t', p = Ps +p, c = c s + c 

The perturbations must satisfy: 

V • u' = 

= -W-^'-p (f3T'-P s C')g 
dV AT „ 2m , 

at H 
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6.7.4 Normalization 



(x,y,z) = H(x*,y*,z*), 
aH 2 



V = 



V* 

IT 



t 



.f, u' = -^u* T = \AT\T* 
H 



C = \AC\C* P ' = ^p* 

k 

The dimensionless equations, are, with * omitted for brevity, 

V-u = 



where 



Ra 
Ra s 
Le 
M 



- Vp - u + (RaT - ^C^j e 3 = 



dT 
~dt 



- w = V 2 T 



V 2 C, 



gf3kHAT 

VK m 

gfckHAC 



UK, 



thermal Rayleigh No. 



= salinity Rayleigh No. 



k 

""rr. 

k s 
n 



Lewis No. 



(6.7.13) 



Note that 
where 



Ra s = RaN Le, 
P S AC 



N 



(3AT 



can be either positive or negative. Taking curl of Eqn. (6.7.15) 

V x u = Ra {{T y - NC y ) e x - (T x - NC X ) e 2 } . 
Taking curl of the above and using V • u = 

V 2 u = -Ra{(T xz -NC xz )e 1 + (T yz -NC yz )e 2 

~ [(T XX + Tyy) ~ N [C XX + Cyy)] 63} . 



The z component is 

d 2 d 2 

V'w = Ra {V\T - NV\C) , V 2 = — + — (6.7.28) 
with the boundary conditions, 

w = T = C = 0, at 2 = 0,1. (6.7.29) 
Consider sinusoidal disturbances 



w \ 




f w 




l-l 


f 


c 




v c 



-iu)t+i(x+imy 



From Eqn. (6.7.28) 

From Eqn. (6.7.16) 

From Eqn. (6.7.17) 

where 



(6.7.30) 



(D 2 - a 2 ) W = -a 2 Ra [T - NCj (6.7.31) 

-iujf -W ={D 2 - a 2 ) f (6.7.32) 

(-iuj{n/a)C - W^j Le = (D 2 - a 2 ) C (6.7.33) 

D = ^- and a 2 = i 2 + m 2 

dz 



The boundary conditions are 

W = f = C = 0, 2 = 0, 1 (6.7.34) 

We can eliminate T and C by applying the operator 

(D 2 -a 2 + iuj) (D 2 -a 2 - iu(n/a)Le) 

on Eqn. (6.7.31) and using Eqns. (6.7.32) and (6.7.33), to get 

(D 2 -a 2 + itu) (D 2 -a 2 + iuj{ri/<j)Le) (D 2 - a 2 ) W 

-a 2 Ra (D 2 - a 2 + iu{n/a)Le) W + a 2 Ra s (D 2 - a 2 + iu) W = (6.7.35) 

The boundary conditions are 

W = (6.7.36) 

D 2 W = see (6.7.31) (6.7.37) 

D 4 W = see (6.7.32) & (6.7.33) (6.7.38) 
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Similarly one can show that D W = on the boundaries. Assume the solution 

W = sinjirz 

we get from Eqn. (6.7.35) the eigenvalue condition, 



a 2 - iuj) (fir 2 + a 2 - iuiAf) 



- (A 2 + u 2 ) (fn 2 + 

+a 2 Ra [fir 2 + a 2 - iujM] - a 2 Ra s [j 2 ir 2 + a 2 - iuj] =0 



where 



M = —Le 
a 



(6.7.39) 



(6.7.40) 



(6.7.41) 



In general u is complex. At the threshold of instability, u = real. From the the real part 
of (6.7.40): 



Ra - Ra s = -u 

a z a z 



(6.7.42) 



From the imaginary part of Eqn. (6.7.40): 

uj [M(aRa) — Ra s ] = u 



° v r J) ' (i + Ao' 

a 2 



(6.7.43) 



Both these equations represent instability thresholds and must hold for real u. 
6.7.5 Monotonic instability 

There is no oscillation and the real part of u is also zero. It follows from (6.7.42) that 

(jV + a 2 f 



Ra — Ra s 



(6.7.44) 



When there is no salinity AC = and only temperature variation in the background, it 
is necessary that AT > (or aRa > , warmer water is at the bottom) for instability. 
Similarly when there is no temperature variation (AT = 0) and only salinity variation, 
AC < ( or Ra s < ) is necessary for instability. With both effects present, we expect 
that Ra — Ra s > corresponds to instability. In particular to find the lowest threshold of 
instability as a function of a, we take j — 1, we take 



d 
da 2 



(Ra - Ra s ) = 



2 2 

a = 7r 



(6.7.45) 



Therefore, the lowest threshold for u = gives 



Ra - Ra, = Air 2 



(6.7.46) 



which is a straight line in Figure 6.7.3. The region of instabilty but lie to the right of this 
line. 
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Figure 6.7.3: Domains of montonic stability and instability 
6.7.6 Oscillatory instability 

We now assume that the real part of lo ^ 0. From Eqn. (6.7.43) 

r 7 - 2 7r 2 + a 2 ) 2 
MRa - Ra s = u Z (1+M) 



The lowest threshold of oscillatory instability is at : j — 1, a 2 — 7T 2 



MRa - Ra s = 4ir 2 (l +J\f) (B) 



(6.7.47) 



If MRa — Ra s > 4ir 2 (l +M), oscillatory instability occurs. The frequency of oscillation u is 
given by Eqn. (6.7.42) 



Mi 



CO 



7T Z 



= 4tt 2 - (Ra - Ra s ) 



Since lo is real, we must have 



47T > Ra - Ra, 



, (A) 



(6.7.48) 



(6.7.49) 



for instability. The equality sign is one boundary. 

The second boundary is defined by Eqn. (6.7.47) which is a straight line that intersects 
with Eqn. (6.7.49) at 



4tt 2 M Ait 2 
( Ra )c = 77 7 and ( Ra s) c 



M- 1 



M-r 



(6.7.50) 



Example: M =2 



{Ra) c = 8tt 2 , (Ra s ) c = An 2 



Hence the wedgedike region (doubly hatched) of oscillatory instability is in the first quadrant 
of Figure 6.7.4. 
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Figure 6.7.4: Domains of stability/instability for nLe = 2. 

Example: TV = 1/2 
The second boundary intersects with the first at 

(Ra) c = -Ati 2 (Ra s ) c = -2 (4vr 2 ) 

In this case the wedgedike region (doubly hatched) of oscillatory instability is in the third 
quadrant of Figure 6.7.5. 



oal 




Figure 6.7.5: Domains of stability/instability for nLe = 1/2. 



6.7.7 Comments on the nonlinear stage 

In a Hele-Shaw experiments Griffiths (1981) showed that the monotonic instability leads to 
fingering. 
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Imhoff & Green have also given laboratory evidence of fingering as a result of such 
instability . The mechanism can be of interest in the transport of 
contaminants in ground water. 

6.7.8 Comments on thermalhaline instabilty in sea water 

Motivated by oceanographic interest, Turner has performed laboratory experiments 
for more rapidly diffusing salt solution on top of denser sucrose solution. 

Unstable growth of the oscillatory instability results in a mixed layer where the salinity 
and temperature are relatively uniform. At the top of this mixed layer, the salinity is 
rather discontinuous because of the weak diffusivity. But the temperature profile has a 
smoother boundary-layer-like transition because of the large diffusivity. The density profile 
is nevertheless continuous. Now this thermal boundary layer is like a Rayleigh -Benard 
(or Horton-Rogers-Lapwood) problem. When the jump become large enough the critical 
Rayleigh number is exceeded so that Benard/Lapwood convection cells appear as sketched 
in Figure 6.7.11. This process continues until another mixed layer develops. Eventually a 
staircase profile is developed, as sketched in figure 6.7.12. 

Field evidences have been found in oceans where there are strong gradients of temperature 
and salinity. It has been observed and recorded in Mediteranian Sea "below the warm 
salty Mediteranian outflow in the eastern North Altlantic and in the western tropical North 
Altlantc" (Schmitt, 1994). This is reflected by the staircase structure of temperature and 
salinity profiles. 
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Figure 6.7.11: Distributions of salinity, temperature, and 
density, generated by beating a iinear saiinity gradient 
from beiow. {Adapted from Kundu, Fiuid Mechanics. 1990). 
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Figure 6.T.1Z: Two types of double-diffusive instabilities, 
{a) Finger instability and (b) Oscillating instability. 
{Adapted from Kundu, Fluid Mechanics, 1990). 



